1 Cluster Analysis

rm(list = ls())

euro <- read.table("/Users/aubrey/Documents/GitHub/MVA_R/europa.txt",header=T, fileEncoding="latin1")

Xeuro <- data.frame(euro[,-1],row.names=euro[,1])
euroscaled <- scale(Xeuro,scale=TRUE,center=TRUE)

dist_eucl2 <- dist(euroscaled, method="euclidean")^2

sl   <- hclust(dist_eucl2, method="single")

com <- hclust(dist_eucl2, method="complete")


# centroid using Mahalanobis dist 
S <- var(euroscaled)
dist_mahal <- matrix(NA,nrow=24,ncol=24)
for(i in 1:24){
  dist_mahal[i,] <- mahalanobis(euroscaled, euroscaled[i,],S)
}
rownames(dist_mahal) <- colnames(dist_mahal) <- labels(dist_eucl2)
dist_mahal <- as.dist(dist_mahal)

cent <- hclust(dist_mahal, method="centroid")

average <- hclust(dist_eucl2, method="average")

#par(mfrow=c(2,2))
plot(sl)

plot(com)

plot(average)

plot(cent)

library(cluster)
#help("kmeans")
set.seed(77)
km4 <- kmeans(x=euroscaled,centers = 4,iter.max = 20,nstart = 10)

#Yields following solution:
km4.cl <- cbind(as.character(euro[order(km4$cluster),1]),sort(km4$cluster))

#Hierarchical Cluster
#Single linkage restricted to arbl and brut with k=4
sl2     <- hclust(dist(euroscaled[,3:4]),method="single")
sl24    <- cutree(sl2,k=4)
sl24.cl <- cbind(as.character(euro[order(sl24),1]),sort(sl24))

# Centroid restricted to arbl and brut with k=4
cent2     <- hclust(dist(euroscaled[,3:4]),method="centroid")
cent24    <- cutree(cent2,k=4)
cent24.cl <- cbind(as.character(euro[order(cent24),1]),sort(cent24))

# Complete Linkage restricted to arbl and brut with k=4 and Mahalanobis dist 
S <- var(euroscaled[,3:4])
dist_mahal <- matrix(NA,nrow=24,ncol=24)
for(i in 1:24){
  dist_mahal[i,] <- mahalanobis(euroscaled[,3:4], euroscaled[i,3:4],S)
}
dist_mahal <- as.dist(dist_mahal)
cl2     <- hclust(dist_mahal, method="complete")
cl24    <- cutree(cl2,k=4)
cl24.cl <- cbind(as.character(euro[order(cl24),1]),sort(cl24))

# kmeans-Clustering restricted to arbl and brut, 10 random initial partitions
set.seed(77)
km2 <- kmeans(x=euroscaled[,3:4],centers=4,iter.max=20,nstart=10)
km2.cl <- cbind(as.character(euro[order(km2$cluster),1]),sort(km2$cluster))

# Graphical comparison
plot_clust <- function(res,main=""){
  plot(euro[res==1,4],euro[res==1,5],type="n",xlim=c(0,37000),ylim=c(0,20),xlab="brut",ylab="arbl",main=main)
  for(i in 1:4){
    text(euro[res==i,4],euro[res==i,5],labels=as.character(euro[res==i,1]),col=i)
  }
}

#par(mfrow=c(2,2))

# kmeans
plot_clust(km2$cluster,main="kmeans")

# single
plot_clust(sl24,main="single")

# centroid
plot_clust(cent24,main="centroid")

# complete
plot_clust(cl24,main="complete")

par(mfrow=c(1,1))

### mclust

library(MASS)
data(geyser)

plot(geyser)

# yielding 3 classes (visually) 


library(mclust)

## Package 'mclust' version 5.4.5
## Type 'citation("mclust")' for citing this R package in publications.

#help(Mclust)
#help(mclustModelNames)


set.seed(0)
m1 <- Mclust(data=geyser, G=3, modelNames = rep("EII", 3))
summary(m1, parameters=TRUE)                                 ## BIC = -3770.529

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust EII (spherical, equal volume) model with 3 components: 
## 
##  log-likelihood   n df       BIC      ICL
##       -1859.609 299  9 -3770.522 -3804.12
## 
## Clustering table:
##   1   2   3 
## 123  77  99 
## 
## Mixing probabilities:
##         1         2         3 
## 0.4097472 0.2644243 0.3258285 
## 
## Means:
##               [,1]      [,2]     [,3]
## waiting  76.425282 87.672625 54.68080
## duration  3.229337  2.620501  4.43386
## 
## Variances:
## [,,1]
##           waiting duration
## waiting  11.24101  0.00000
## duration  0.00000 11.24101
## [,,2]
##           waiting duration
## waiting  11.24101  0.00000
## duration  0.00000 11.24101
## [,,3]
##           waiting duration
## waiting  11.24101  0.00000
## duration  0.00000 11.24101

plot(m1, what="classification")

plot(m1, what="uncertainty")

plot(m1, what="density")               #

second model

m2 <- Mclust(data=geyser, G=3, modelNames = rep("VVI", 3))
summary(m2, parameters=TRUE)                                 ## BIC = -2813.52

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VVI (diagonal, varying volume and shape) model with 3 components: 
## 
##  log-likelihood   n df       BIC       ICL
##       -1369.343 299 14 -2818.493 -2831.405
## 
## Clustering table:
##   1   2   3 
##  95 107  97 
## 
## Mixing probabilities:
##         1         2         3 
## 0.3117014 0.3549744 0.3333242 
## 
## Means:
##               [,1]      [,2]      [,3]
## waiting  78.347255 83.257569 55.018887
## duration  4.084144  1.991457  4.442714
## 
## Variances:
## [,,1]
##           waiting  duration
## waiting  40.92059 0.0000000
## duration  0.00000 0.1033852
## [,,2]
##           waiting   duration
## waiting  43.76243 0.00000000
## duration  0.00000 0.08715797
## [,,3]
##           waiting  duration
## waiting  31.33197 0.0000000
## duration  0.00000 0.1237322

plot(m2, what="classification")

plot(m2, what="uncertainty")

plot(m2, what="density")                                     ## improvement subject to m1

third model

m3 <- Mclust(data=geyser, G=3, modelNames = rep("VVV", 3))
summary(m3, parameters=TRUE)                                 ## BIC = -3059.029

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VVV (ellipsoidal, varying volume, shape, and orientation) model
## with 3 components: 
## 
##  log-likelihood   n df       BIC       ICL
##       -1364.937 299 17 -2826.782 -2840.339
## 
## Clustering table:
##   1   2   3 
## 101 101  97 
## 
## Mixing probabilities:
##         1         2         3 
## 0.3409570 0.3376228 0.3214202 
## 
## Means:
##               [,1]      [,2]      [,3]
## waiting  78.271156 83.129899 54.634820
## duration  4.043017  1.945981  4.434418
## 
## Variances:
## [,,1]
##             waiting   duration
## waiting  47.6582109 -0.7605575
## duration -0.7605575  0.1819945
## [,,2]
##            waiting    duration
## waiting  44.399235 -0.27259302
## duration -0.272593  0.04874916
## [,,3]
##              waiting    duration
## waiting  27.99788820 -0.03423172
## duration -0.03423172  0.12372777

plot(m3, what="classification")                              ## does not fit well

plot(m3, what="uncertainty")

plot(m3, what="density")                                     ## does not fit well

## estimated posteriori over all obs:
round(m2$z, 3)

##      [,1]  [,2]  [,3]
## 1   1.000 0.000 0.000
## 2   0.000 1.000 0.000
## 3   0.008 0.000 0.992
## 4   1.000 0.000 0.000
## 5   0.999 0.000 0.001
## 6   0.000 1.000 0.000
## 7   0.014 0.000 0.986
## 8   1.000 0.000 0.000
## 9   0.000 1.000 0.000
## 10  0.000 0.000 1.000
## 11  0.000 1.000 0.000
## 12  0.000 0.000 1.000
## 13  0.000 1.000 0.000
## 14  0.000 0.000 1.000
## 15  1.000 0.000 0.000
## 16  0.000 1.000 0.000
## 17  0.005 0.000 0.995
## 18  0.000 1.000 0.000
## 19  0.082 0.000 0.918
## 20  1.000 0.000 0.000
## 21  0.000 1.000 0.000
## 22  0.000 0.000 1.000
## 23  0.000 1.000 0.000
## 24  0.014 0.000 0.986
## 25  0.000 1.000 0.000
## 26  0.008 0.000 0.992
## 27  0.009 0.991 0.000
## 28  0.609 0.000 0.391
## 29  0.999 0.000 0.001
## 30  1.000 0.000 0.000
## 31  0.999 0.000 0.001
## 32  1.000 0.000 0.000
## 33  0.000 1.000 0.000
## 34  0.046 0.000 0.954
## 35  0.000 1.000 0.000
## 36  0.009 0.000 0.991
## 37  0.000 1.000 0.000
## 38  0.000 0.000 1.000
## 39  0.000 1.000 0.000
## 40  0.000 0.000 1.000
## 41  0.000 1.000 0.000
## 42  0.001 0.000 0.999
## 43  0.000 1.000 0.000
## 44  0.001 0.000 0.999
## 45  0.000 1.000 0.000
## 46  0.000 0.000 1.000
## 47  0.000 1.000 0.000
## 48  0.021 0.000 0.979
## 49  0.000 1.000 0.000
## 50  0.000 0.000 1.000
## 51  0.000 1.000 0.000
## 52  0.044 0.000 0.956
## 53  0.000 1.000 0.000
## 54  0.000 0.000 1.000
## 55  0.000 1.000 0.000
## 56  0.010 0.000 0.990
## 57  1.000 0.000 0.000
## 58  0.984 0.000 0.015
## 59  0.999 0.000 0.001
## 60  0.983 0.000 0.017
## 61  0.000 1.000 0.000
## 62  0.000 0.000 1.000
## 63  0.000 1.000 0.000
## 64  0.008 0.000 0.992
## 65  0.000 1.000 0.000
## 66  0.079 0.000 0.921
## 67  0.000 1.000 0.000
## 68  0.000 0.000 1.000
## 69  1.000 0.000 0.000
## 70  0.000 1.000 0.000
## 71  0.032 0.000 0.968
## 72  0.000 1.000 0.000
## 73  0.001 0.000 0.999
## 74  1.000 0.000 0.000
## 75  0.988 0.000 0.012
## 76  0.000 1.000 0.000
## 77  0.018 0.000 0.982
## 78  1.000 0.000 0.000
## 79  0.981 0.000 0.019
## 80  1.000 0.000 0.000
## 81  1.000 0.000 0.000
## 82  0.000 1.000 0.000
## 83  0.025 0.000 0.975
## 84  0.458 0.542 0.000
## 85  0.670 0.000 0.330
## 86  0.000 1.000 0.000
## 87  0.000 0.000 1.000
## 88  0.000 1.000 0.000
## 89  0.000 0.000 1.000
## 90  0.000 1.000 0.000
## 91  0.000 0.000 1.000
## 92  0.000 1.000 0.000
## 93  0.004 0.000 0.996
## 94  0.000 1.000 0.000
## 95  0.000 0.000 1.000
## 96  0.000 1.000 0.000
## 97  0.000 0.000 1.000
## 98  0.000 1.000 0.000
## 99  0.004 0.000 0.996
## 100 0.000 1.000 0.000
## 101 0.000 0.000 1.000
## 102 0.000 1.000 0.000
## 103 0.000 0.000 1.000
## 104 0.000 1.000 0.000
## 105 0.014 0.000 0.986
## 106 0.000 1.000 0.000
## 107 0.000 0.000 1.000
## 108 0.000 1.000 0.000
## 109 0.000 0.000 1.000
## 110 0.215 0.785 0.000
## 111 0.557 0.000 0.443
## 112 0.999 0.000 0.001
## 113 0.000 1.000 0.000
## 114 0.001 0.000 0.999
## 115 0.000 1.000 0.000
## 116 0.000 0.000 1.000
## 117 0.000 1.000 0.000
## 118 0.000 0.000 1.000
## 119 0.000 1.000 0.000
## 120 0.006 0.000 0.994
## 121 0.000 1.000 0.000
## 122 0.001 0.000 0.999
## 123 0.996 0.000 0.004
## 124 1.000 0.000 0.000
## 125 0.000 1.000 0.000
## 126 0.025 0.000 0.975
## 127 1.000 0.000 0.000
## 128 1.000 0.000 0.000
## 129 0.975 0.000 0.025
## 130 1.000 0.000 0.000
## 131 1.000 0.000 0.000
## 132 1.000 0.000 0.000
## 133 0.000 1.000 0.000
## 134 0.000 0.000 1.000
## 135 1.000 0.000 0.000
## 136 1.000 0.000 0.000
## 137 0.994 0.000 0.006
## 138 0.998 0.000 0.002
## 139 0.000 1.000 0.000
## 140 0.033 0.000 0.967
## 141 1.000 0.000 0.000
## 142 0.946 0.000 0.054
## 143 1.000 0.000 0.000
## 144 1.000 0.000 0.000
## 145 1.000 0.000 0.000
## 146 0.999 0.000 0.001
## 147 0.000 1.000 0.000
## 148 0.074 0.000 0.926
## 149 0.000 1.000 0.000
## 150 0.000 0.000 1.000
## 151 0.000 1.000 0.000
## 152 0.000 0.000 1.000
## 153 0.000 1.000 0.000
## 154 0.000 0.000 1.000
## 155 0.000 1.000 0.000
## 156 0.000 0.000 1.000
## 157 0.000 1.000 0.000
## 158 1.000 0.000 0.000
## 159 1.000 0.000 0.000
## 160 0.991 0.000 0.009
## 161 0.984 0.000 0.016
## 162 1.000 0.000 0.000
## 163 1.000 0.000 0.000
## 164 0.000 1.000 0.000
## 165 0.002 0.000 0.998
## 166 0.000 1.000 0.000
## 167 0.000 0.000 1.000
## 168 0.000 1.000 0.000
## 169 0.016 0.000 0.984
## 170 0.031 0.969 0.000
## 171 0.000 0.000 1.000
## 172 1.000 0.000 0.000
## 173 0.986 0.000 0.014
## 174 0.000 1.000 0.000
## 175 0.002 0.000 0.998
## 176 0.000 1.000 0.000
## 177 0.001 0.000 0.999
## 178 0.000 1.000 0.000
## 179 0.052 0.000 0.948
## 180 1.000 0.000 0.000
## 181 0.000 1.000 0.000
## 182 0.000 0.000 1.000
## 183 0.000 1.000 0.000
## 184 0.000 0.000 1.000
## 185 1.000 0.000 0.000
## 186 1.000 0.000 0.000
## 187 0.601 0.000 0.399
## 188 0.000 1.000 0.000
## 189 0.000 0.000 1.000
## 190 0.000 1.000 0.000
## 191 0.011 0.000 0.989
## 192 0.000 1.000 0.000
## 193 0.000 0.000 1.000
## 194 0.000 1.000 0.000
## 195 0.000 0.000 1.000
## 196 0.994 0.006 0.000
## 197 0.996 0.000 0.004
## 198 0.000 1.000 0.000
## 199 0.031 0.000 0.969
## 200 0.000 1.000 0.000
## 201 0.000 0.000 1.000
## 202 0.000 1.000 0.000
## 203 0.000 0.000 1.000
## 204 1.000 0.000 0.000
## 205 0.000 1.000 0.000
## 206 0.515 0.000 0.485
## 207 0.999 0.000 0.001
## 208 0.000 1.000 0.000
## 209 0.803 0.000 0.197
## 210 1.000 0.000 0.000
## 211 0.712 0.000 0.288
## 212 0.000 1.000 0.000
## 213 0.020 0.000 0.980
## 214 0.000 1.000 0.000
## 215 0.001 0.000 0.999
## 216 0.000 1.000 0.000
## 217 0.000 0.000 1.000
## 218 0.000 1.000 0.000
## 219 0.984 0.000 0.016
## 220 0.992 0.000 0.008
## 221 0.000 1.000 0.000
## 222 0.000 0.000 1.000
## 223 1.000 0.000 0.000
## 224 1.000 0.000 0.000
## 225 0.998 0.000 0.002
## 226 0.984 0.000 0.016
## 227 1.000 0.000 0.000
## 228 0.997 0.000 0.003
## 229 0.000 1.000 0.000
## 230 0.000 0.000 1.000
## 231 0.000 1.000 0.000
## 232 0.000 0.000 1.000
## 233 0.000 1.000 0.000
## 234 0.000 0.000 1.000
## 235 1.000 0.000 0.000
## 236 1.000 0.000 0.000
## 237 1.000 0.000 0.000
## 238 0.000 1.000 0.000
## 239 0.006 0.000 0.994
## 240 1.000 0.000 0.000
## 241 0.000 1.000 0.000
## 242 0.000 0.000 1.000
## 243 0.157 0.843 0.000
## 244 0.991 0.000 0.009
## 245 0.000 1.000 0.000
## 246 0.001 0.000 0.999
## 247 0.998 0.002 0.000
## 248 0.000 1.000 0.000
## 249 0.124 0.000 0.876
## 250 0.000 1.000 0.000
## 251 0.001 0.000 0.999
## 252 1.000 0.000 0.000
## 253 0.982 0.000 0.018
## 254 0.000 1.000 0.000
## 255 0.000 0.000 1.000
## 256 0.000 1.000 0.000
## 257 0.003 0.000 0.997
## 258 1.000 0.000 0.000
## 259 0.580 0.000 0.420
## 260 1.000 0.000 0.000
## 261 0.992 0.000 0.008
## 262 0.000 1.000 0.000
## 263 0.001 0.000 0.999
## 264 1.000 0.000 0.000
## 265 1.000 0.000 0.000
## 266 0.000 1.000 0.000
## 267 0.000 0.000 1.000
## 268 0.000 1.000 0.000
## 269 0.000 0.000 1.000
## 270 0.001 0.999 0.000
## 271 0.025 0.000 0.975
## 272 1.000 0.000 0.000
## 273 0.000 1.000 0.000
## 274 0.005 0.000 0.995
## 275 0.000 1.000 0.000
## 276 0.000 0.000 1.000
## 277 1.000 0.000 0.000
## 278 1.000 0.000 0.000
## 279 1.000 0.000 0.000
## 280 1.000 0.000 0.000
## 281 1.000 0.000 0.000
## 282 0.999 0.000 0.001
## 283 0.998 0.000 0.002
## 284 0.000 1.000 0.000
## 285 0.670 0.000 0.330
## 286 0.000 1.000 0.000
## 287 0.000 0.000 1.000
## 288 0.000 1.000 0.000
## 289 0.000 0.000 1.000
## 290 0.000 1.000 0.000
## 291 0.000 0.000 1.000
## 292 0.000 1.000 0.000
## 293 0.000 0.000 1.000
## 294 0.000 1.000 0.000
## 295 0.000 0.000 1.000
## 296 0.000 1.000 0.000
## 297 0.014 0.000 0.986
## 298 1.000 0.000 0.000
## 299 0.000 1.000 0.000

hist(round(m2$z, 3))

## classification:
m2$classification

##   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18 
##   1   2   3   1   1   2   3   1   2   3   2   3   2   3   1   2   3   2 
##  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36 
##   3   1   2   3   2   3   2   3   2   1   1   1   1   1   2   3   2   3 
##  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54 
##   2   3   2   3   2   3   2   3   2   3   2   3   2   3   2   3   2   3 
##  55  56  57  58  59  60  61  62  63  64  65  66  67  68  69  70  71  72 
##   2   3   1   1   1   1   2   3   2   3   2   3   2   3   1   2   3   2 
##  73  74  75  76  77  78  79  80  81  82  83  84  85  86  87  88  89  90 
##   3   1   1   2   3   1   1   1   1   2   3   2   1   2   3   2   3   2 
##  91  92  93  94  95  96  97  98  99 100 101 102 103 104 105 106 107 108 
##   3   2   3   2   3   2   3   2   3   2   3   2   3   2   3   2   3   2 
## 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 
##   3   2   1   1   2   3   2   3   2   3   2   3   2   3   1   1   2   3 
## 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 
##   1   1   1   1   1   1   2   3   1   1   1   1   2   3   1   1   1   1 
## 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 
##   1   1   2   3   2   3   2   3   2   3   2   3   2   1   1   1   1   1 
## 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 
##   1   2   3   2   3   2   3   2   3   1   1   2   3   2   3   2   3   1 
## 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 
##   2   3   2   3   1   1   1   2   3   2   3   2   3   2   3   1   1   2 
## 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 
##   3   2   3   2   3   1   2   1   1   2   1   1   1   2   3   2   3   2 
## 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 
##   3   2   1   1   2   3   1   1   1   1   1   1   2   3   2   3   2   3 
## 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 
##   1   1   1   2   3   1   2   3   2   1   2   3   1   2   3   2   3   1 
## 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 
##   1   2   3   2   3   1   1   1   1   2   3   1   1   2   3   2   3   2 
## 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 
##   3   1   2   3   2   3   1   1   1   1   1   1   1   2   1   2   3   2 
## 289 290 291 292 293 294 295 296 297 298 299 
##   3   2   3   2   3   2   3   2   3   1   2

## Uncertainty subject to partition:
round(m2$uncertainty, 3)

##     1     2     3     4     5     6     7     8     9    10    11    12 
## 0.000 0.000 0.008 0.000 0.001 0.000 0.014 0.000 0.000 0.000 0.000 0.000 
##    13    14    15    16    17    18    19    20    21    22    23    24 
## 0.000 0.000 0.000 0.000 0.005 0.000 0.082 0.000 0.000 0.000 0.000 0.014 
##    25    26    27    28    29    30    31    32    33    34    35    36 
## 0.000 0.008 0.009 0.391 0.001 0.000 0.001 0.000 0.000 0.046 0.000 0.009 
##    37    38    39    40    41    42    43    44    45    46    47    48 
## 0.000 0.000 0.000 0.000 0.000 0.001 0.000 0.001 0.000 0.000 0.000 0.021 
##    49    50    51    52    53    54    55    56    57    58    59    60 
## 0.000 0.000 0.000 0.044 0.000 0.000 0.000 0.010 0.000 0.016 0.001 0.017 
##    61    62    63    64    65    66    67    68    69    70    71    72 
## 0.000 0.000 0.000 0.008 0.000 0.079 0.000 0.000 0.000 0.000 0.032 0.000 
##    73    74    75    76    77    78    79    80    81    82    83    84 
## 0.001 0.000 0.012 0.000 0.018 0.000 0.019 0.000 0.000 0.000 0.025 0.458 
##    85    86    87    88    89    90    91    92    93    94    95    96 
## 0.330 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.004 0.000 0.000 0.000 
##    97    98    99   100   101   102   103   104   105   106   107   108 
## 0.000 0.000 0.004 0.000 0.000 0.000 0.000 0.000 0.014 0.000 0.000 0.000 
##   109   110   111   112   113   114   115   116   117   118   119   120 
## 0.000 0.215 0.443 0.001 0.000 0.001 0.000 0.000 0.000 0.000 0.000 0.006 
##   121   122   123   124   125   126   127   128   129   130   131   132 
## 0.000 0.001 0.004 0.000 0.000 0.025 0.000 0.000 0.025 0.000 0.000 0.000 
##   133   134   135   136   137   138   139   140   141   142   143   144 
## 0.000 0.000 0.000 0.000 0.006 0.002 0.000 0.033 0.000 0.054 0.000 0.000 
##   145   146   147   148   149   150   151   152   153   154   155   156 
## 0.000 0.001 0.000 0.074 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 
##   157   158   159   160   161   162   163   164   165   166   167   168 
## 0.000 0.000 0.000 0.009 0.016 0.000 0.000 0.000 0.002 0.000 0.000 0.000 
##   169   170   171   172   173   174   175   176   177   178   179   180 
## 0.016 0.031 0.000 0.000 0.014 0.000 0.002 0.000 0.001 0.000 0.052 0.000 
##   181   182   183   184   185   186   187   188   189   190   191   192 
## 0.000 0.000 0.000 0.000 0.000 0.000 0.399 0.000 0.000 0.000 0.011 0.000 
##   193   194   195   196   197   198   199   200   201   202   203   204 
## 0.000 0.000 0.000 0.006 0.004 0.000 0.031 0.000 0.000 0.000 0.000 0.000 
##   205   206   207   208   209   210   211   212   213   214   215   216 
## 0.000 0.485 0.001 0.000 0.197 0.000 0.288 0.000 0.020 0.000 0.001 0.000 
##   217   218   219   220   221   222   223   224   225   226   227   228 
## 0.000 0.000 0.016 0.008 0.000 0.000 0.000 0.000 0.002 0.016 0.000 0.003 
##   229   230   231   232   233   234   235   236   237   238   239   240 
## 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.006 0.000 
##   241   242   243   244   245   246   247   248   249   250   251   252 
## 0.000 0.000 0.157 0.009 0.000 0.001 0.002 0.000 0.124 0.000 0.001 0.000 
##   253   254   255   256   257   258   259   260   261   262   263   264 
## 0.018 0.000 0.000 0.000 0.003 0.000 0.420 0.000 0.008 0.000 0.001 0.000 
##   265   266   267   268   269   270   271   272   273   274   275   276 
## 0.000 0.000 0.000 0.000 0.000 0.001 0.025 0.000 0.000 0.005 0.000 0.000 
##   277   278   279   280   281   282   283   284   285   286   287   288 
## 0.000 0.000 0.000 0.000 0.000 0.001 0.002 0.000 0.330 0.000 0.000 0.000 
##   289   290   291   292   293   294   295   296   297   298   299 
## 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.014 0.000 0.000

table(round(m2$uncertainty, 3))

## 
##     0 0.001 0.002 0.003 0.004 0.005 0.006 0.008 0.009  0.01 0.011 0.012 
##   206    19     6     2     4     2     4     5     4     1     1     1 
## 0.014 0.016 0.017 0.018 0.019  0.02 0.021 0.025 0.031 0.032 0.033 0.044 
##     5     5     1     2     1     1     1     4     2     1     1     1 
## 0.046 0.052 0.054 0.074 0.079 0.082 0.124 0.157 0.197 0.215 0.288  0.33 
##     1     1     1     1     1     1     1     1     1     1     1     2 
## 0.391 0.399  0.42 0.443 0.458 0.485 
##     1     1     1     1     1     1

four components

Mclust(data=geyser)

## 'Mclust' model object: (VVI,4) 
## 
## Available components: 
##  [1] "call"           "data"           "modelName"      "n"             
##  [5] "d"              "G"              "BIC"            "bic"           
##  [9] "loglik"         "df"             "hypvol"         "parameters"    
## [13] "z"              "classification" "uncertainty"

m4 <- Mclust(data=geyser, G=4, modelNames = rep("VVI", 4))
summary(m4, parameters=TRUE)                                 ## BIC = -2768.548

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VVI (diagonal, varying volume and shape) model with 4 components: 
## 
##  log-likelihood   n df       BIC       ICL
##        -1330.13 299 19 -2768.568 -2798.746
## 
## Clustering table:
##  1  2  3  4 
## 90 17 98 94 
## 
## Mixing probabilities:
##          1          2          3          4 
## 0.29283469 0.07968364 0.33392296 0.29355871 
## 
## Means:
##               [,1]      [,2]      [,3]      [,4]
## waiting  78.245444 81.693185 55.045717 83.495263
## duration  4.124353  2.585291  4.443712  1.918517
## 
## Variances:
## [,,1]
##           waiting   duration
## waiting  41.31804 0.00000000
## duration  0.00000 0.07151279
## [,,2]
##          waiting duration
## waiting  33.8808  0.00000
## duration  0.0000  0.44011
## [,,3]
##           waiting  duration
## waiting  31.65528 0.0000000
## duration  0.00000 0.1239616
## [,,4]
##           waiting   duration
## waiting  45.44063 0.00000000
## duration  0.00000 0.01695964

plot(m4, what="classification")

plot(m4, what="uncertainty")

plot(m4, what="density")

Block_B1

Dingyi Lai

1/8/2021

1 Cluster Analysis